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1.  Introduction 


In  texture  segmentation  it  is  appropriate  to  process  the 
image  with  operators  of  many  sizes.  This  is  computationally 
expensive  because  many  computational  steps  are  required  for 
each  location  in  the  image. 

In  [1]  ,  Burt  describes  an  elegant  hierarchical  discrete 
correlation  (HDC)  method  for  computing  the  correlation  between 
an  image  and  a  kernel.  In  this  approach  correlations  at  upper 
levels  can  be  computed  as  the  weighted  sums  of  correlations 
with  narrower  kernels  at  lower  levels.  Kernels  which  can  be 
computed  hierarchically  in  this  way  closely  approximate  the 
Gaussian  probability  distribution.  This  means  that  correla¬ 
tion  is  equivalent  to  low-pass  filtering. 

The  principal  advantage  of  the  HDC  method  is  that  it  is 
computationally  more  efficient  than  the  direct  correlation 
and  FFT  methods.  In  addition,  correlations  for  a  set  of 
scaled  kernels  are  computed  at  once,  without  any  need  to 
construct  and  store  large  kernels  or  kernels  of  different 
shapes  and  sizes.  Samples  of  the  correlations  obtained  at 
nearby  image  positions  can  be  summed  to  obtain  e.g.  band-pass 
Laplacian  ("Mexican  hat")  operators.  Mexican  hat  filtered 
images  are  essentially  the  same  as  the  differences  of 
Gaussian  filtered  channels  observed  in  the  human  visual 
system  [2]. 

Band-pass  filtering  responds  to  details  of  the  image  which 
contain  a  limited  range  of  spatial  frequencies.  This  reflects 
the  coarseness  of  the  texture.  Thus,  features  derived  from 
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band-pass  filtered  images  on  different  frequency  bands 
might  be  useful  for  texture  classification.  The  local 
computation  of  these  textural  properties  could  be  an  aid 
to  image  segmentation  using  cooperative  hierarchical  com¬ 
putation  [3] . 

In  this  note  we  apply  Mexican  hat  filtering  to  180 
LANDSAT  imagery  samples  belonging  to  three  geological 
terrain  types.  The  averages  of  the  squared  outputs  computed 
over  band-pass  filtered  images  are  used  as  texture  features 
in  a  Fisher  linear  discriminant  classification.  The  results 
are  compared  to  the  results  obtained  by  Weszka  et  al.  [4]. 

Nondirectional  Mexican  hat  features  did  not  yield  as 
good  classification  results  as  some  of  the  statistical  tex¬ 
ture  features  using  directional  information.  Their  performance 
was  more  comparable  to  that  of  features  derived  from  the 
Fourier  power  spectrum.  In  any  case,  Mexican  hat  features 
seem  to  give  a  good  measure  of  texture  coarseness.  Combined 
with  other  measures,  they  should  be  useful  for  texture  segmen¬ 


tation. 


2.1  Mexican  hat  filtering 

The  band-limited  Laplacian  ("Mexican  hat")  may  be  formally 
defined  as  the  Laplacian  of  a  Gaussian  A2G,  but  the  result 
is  generally  approximated  as  the  difference  between  two 
Gaussian  functions  which  have  different  standard  deviations 
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Typically,  the  ratio  ^ ai  in  ran9e  1,5  to  2.5. 

Burt  defined  three  types  of  HDC.  In  type  1  an  odd  number 
and  in  type  2  an  even  number  of  correlations  with  small 
kernels  are  summed  to  obtain  the  correlation  with  the  next 
larger  kernel,  while  type  3  is  for  fractional  values.  He  also 
defined  a  reduced  form  of  HDC  in  which  the  number  of  samples 
is  reduced  by  a  factor  of  r2  from  level  to  level. 

Let  g^(x,y)  be  the  correlation  function  at  hierarchical 
level  £.  It  was  obtained  from  the  original  image  f(x,y)  through 
l  recursions  of  a  correlation-like  operation  using  the  weighting 
function  w(x,y) .  The  sample  distance  grows  geometrically  by 
the  factor  r  from  level  to  level;  thus  r  is  said  to  be  the 


order  of  the  HDC. 
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A  band-limited  Laplacian  with  a  ratio  o2/o  “ 
be  obtained  from  the  type  1  HDC  of  order  r  by  subtracting 

g£+l(x'y)  from  g£(x,y)*  If  the  reduced  form  of  the  HDC 
is  used,  samples  of  9£+1(x'y)  are  not  computed  for  every 

sample  of  g„ (x,y) .  Missing  samples  are  obtained  by  apply- 

ing  w (x,y)  to  the  neighborhood  of  g^(x,y). 


To  obtain  ratios  a2/o1<  r  two  HDC's  (of  any  type)  are 
computed  ’ith  different  generating  kernels.  The  Laplacian 
operator  is  then  a  difference  between  corresponding  samples 
in  the  two  HDC  *  s . 

For  r'  f^rst  an  HDC  of  any  type  is  obtained 

to  form  the  central  Gaussian.  Then  the  surround  Gaussian 
is  obtained  by  applying  the  generating  kernel  in  reverse. 
For  type  1  we  have 
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The  sums  in  this  expression  are  understood  to  include  only 

n  -i  n  -j 

those  terms  for  which  ~ —  and  are  integer  valued.  The 

equal  contribution  constraint  on  w  ensures  that  the  sum  will 
have  a  total  weight  at  1/r2,  hence  the  r2  factor  normalizes 
the  sum  in  the  above  definition. 

A  one-dimensional  example  using  this  procedure  is  shown 
in  Figure  1  [1].  Here  a  ratio  o2/a1  -  2.5  is  obtained  with 


a  type  1  HDC,  order  r=2,  kernel  width  k»5,  and  with  a  separ¬ 
able  generating  kernel  with  weight  vector  w  *  (.05, .25, .4, .25, .05) 
The  Mexican  hat  filtering  algorithm  used  in  our  study  was 
based  on  this  third  procedure. 


2.2  Features  and  classifier 

As  a  result  of  the  hierarchical  computation,  we  have  an 
order  N  pyramid,  in  which  the  lowest  level  represents  the 
highest  spatial  frequency  band.  Its  dimensions  are  the  same 
as  in  the  original  image.  Image  array  dimensions  are  de¬ 
creased  by  half  from  level  to  level,  and  these  levels  repre¬ 
sent  lower  frequency  bands. 

Let  L.  .  ,  be  the  value  of  the  filtered  image  at  point 
(i,j)  on  level  l.  Texture  feature  LAP  is  defined  as  an 
average  of  the  squared  values  computed  over  the  filtered  image 

^  2exp(N+l-Jl)  2exp(N+l-Jt) 

LAP(£)  =  2exp (N+1-& j  2exp(N+l-£)  ?=1  Li,j,£ 

for  i«l, 2,..., N.  For  a  64  by  64  image,  for  example,  we  have 
six  features.  The  first  of  these  is  computed  over  a  64  by  64 
image  and  the  last  one  over  a  2  by  2  image. 

The  classifier  used  in  this  study  was  the  same  Fisher 
linear  discriminant  classifier  as  in  [4]. 
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2 . 3  Image  data 


The  same  imagery  as  in  the  main  study  of  [4]  was  used. 

It  comprises  a  set  of  180  LANDSAT  terrain  samples,  belonging 
to  three  geological  terrain  types:  Mississippian  limestone 
and  shale;  Lower  Pennsylvanian  shale,  and  Pennsylvanian  sand¬ 
stone  and  shale  (labeled  A,  B  and  C) .  The  gray  scale  has 
been  modified  to  cover  just  64  gray  levels  and  histogram 
flattening  has  been  performed  on  each  of  these  180  terrain 
samples  to  remove  effects  of  unequal  brightness  and  contrast. 


3.  Experiments 

First  a  pilot  study  was  performed  by  using  different 
generating  kernels  in  filtering.  The  generating  kernels  used 
were: 

wl  =  (.05, .25, .4, .25, .05) ,  w2  =  ( . 13 , . 37 , . 37 , . 13)  and 
w3  =  (.2,. 6,. 2).  Here  wl  and  w2  are  the  most  Gaussian-like 
kernels  for  kernel  widths  5  and  4.  The  w3  kernel  is  not 
Gaussian-like,  but  it  was  used  in  this  study  because  it  accen¬ 
tuates  high  spatial  frequencies.  Using  these  kernels  first 
ten  terrain  samples  of  each  type  were  processed.  The  results 
were  quite  similar  for  each  of  these  kernels.  It  seemed  that 
none  of  these  was  able  to  discriminate  very  well  terrain  types 
A  and  B,  while  each  of  them  was  able  to  discriminate  type  C. 
Kernel  w2  was  selected  for  the  main  study  because  it  has 
smaller  border  effects  than  wl  at  upper  pyramid  levels. 

In  the  main  study  all  of  the  180  terrain  samples  were 
processed  using  kernel  w2.  As  a  result  we  got  six  feature 
values  for  each  terrain  sample,  each  of  them  corresponding  to 
a  different  spatial  frequency  band.  The  distributions  for  some 
of  these  features  are  presented  in  Figure  2. 

The  Fisher  linear  discriminant  classification  method  was 
used  to  classify  the  terrain  samples  into  three  classes  by 
using  single  features  and  pairs  of  features.  The  sixth  fea¬ 
ture  was  not  used  in  classification  because  its  values  seemed 
to  be  randomly  distributed,  having  no  correlation  with  the 
different  terrain  types.  Classification  results  for  single 
features  are  presented  in  Table  la  and  for  pairs  of  features 


in  Table  lb. 


.  Discussion  and  conclusions 

A  direct  comparison  to  results  of  Weszka  et  al.  [4]  is 
difficult  because  our  features  have  no  directional  component. 

In  [4] ,  the  features  were  computed  for  all  pixel  sizes  (fre¬ 
quency  bands)  in  four  directions,  and  all  features  and  feature 
pair  combinations  were  used  in  classification.  Neither  are 
the  frequency  bands  used  in  these  studies  exactly  equal. 
However,  by  computing  the  averages  of  the  classification  re¬ 
sults  of  the  earlier  study  for  each  frequency  band,  we  can 
have  a  rough  estimate  of  the  performance  of  the  LAP  features. 

In  Tables  2a  and  2b  are  presented  the  averages  of  the  classi¬ 
fication  results  in  [4]  for  the  Fourier  and  "CONTRAST"  (com¬ 
puted  on  single  points  or  on  averages)  features.  The  best 
classifications  for  each  size  are  presented  in  brackets. 

The  following  can  be  seen  from  these  results: 

1.  Features  based  on  Mexican  hat  (LAP)  filtering  do  not 
quite  as  well  as  features  based  on  second-order  statis¬ 
tics,  and  do  about  equally  well  as  Fourier  features. 

The  absence  of  directional  information  reduces  the 
performance  of  the  LAP  features.  However,  their  per¬ 
formance  does  not  degrade  as  rapidly  with  increasing 
size  as  does  that  of  the  Fourier  and  statistical 
features . 

2.  Results  for  terrain  type  C  are  very  good.  Using  a 
single  feature  only  one  of  the  type  C  samples  was  mis- 
classified.  This  result  is  better  than  the  results 

of  most  of  the  experiments  using  feature  pairs  in  [4j. 


3. 


Discrimination  of  terrain  types  A  and  B  is  quite 
poor.  The  distributions  of  features  at  level  2 
(Figure  2)  show  that  the  feature  values  of  the 
misclassified  samples  are  quite  different  from  the 
values  of  the  other  samples  belonging  to  the  same 
terrain  type. 

These  results  indicate  that  the  LAP  features  give  a  good  mea¬ 
sure  of  texture  coarseness,  but  in  order  to  get  better  classi¬ 
fication  accuracy  some  additional  measures  are  needed.  Prob¬ 
ably  these  samples  are  more  appropriately  modeled  statistically 
in  the  space  domain,  rather  than  as  sums  of  sinusoids,  and 
statistical  features  captured  the  essential  differences  more 
effectively  [4].  This  is  also  in  agreement  with  the  results 
of  Conners  and  Harlow  [5] . 
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Misclassified 

Feature 

number 

Correctly 

classified 

A's 

B's 

C's 

1 

121 

33 

25 

1 

2 

85 

16 

48 

31 

3 

122 

30 

23 

5 

4 
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38 

16 

21 

5 

80 

48 

19 
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Table  la.  Results  using  single  features 

Misclassified 


Feature 

pair 

Correctly 

classified 

A'S 

B’s 

C's 

1,2 

147 

15 

16 

2 

1,3 
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27 

26 

1 

1,4 

137 

21 

20 

2 

1.5 
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24 

22 

2 

2,3 

140 

17 

17 

6 

2,4 

131 

20 

21 

8 

2,5 

110 

18 

30 

22 

3,4 

131 

25 

19 

5 

3,5 

126 

28 

21 

5 

4,5 

107 

38 

15 

20 

Table  lb. 

Results  using 

pairs  of 

features 

CONTRAST 


CONTRAST  COOC  (average 

Feature 

pair 

Fourier 

COOC  (points) 

gray  levels) 

LAP 

Size 

1, 

Size 

1 

142.8  (158) 

151.5(166) 

151.5  (166) 

- 

Size 

1, 

Size 

2 

139.3(151) 

146.3(166) 

143.5(158) 

147 

Size 

1, 

Size 

3 

133.1(144) 

138.6  (153) 

139.1(155) 

126 

Size 

1, 

Size 

4 

128.4 (137) 

130.7  (140) 

138.6 (152) 

137 

Size 

2, 

Size 

2 

116.5  (135) 

141.0  (154) 

118.8  (142) 

- 

Size 

2, 

Size 

3 

117.4 (134) 

132.1  (148) 

130.6  (151) 

140 

Size 

2, 

Size 

4 

107.1(126) 

126.3  (140) 

129.5  (146) 

131 

Size 

3, 

Size 

3 

117.2 (122) 

119.5 (133) 

129.8  (144) 

- 

Size 

3, 

Size 

4 

109.7(120) 

109.6  (126) 

128.8  (139) 

131 

Size 

4, 

Size 

4 

98.7(109) 

93.7(110) 

121.5  (130) 

- 

Table 

•  2a.  Average  values  of 

classification 

results  for 

pairs 

of 

features . 

CONTRAST 

CONTRAST 

Feature 

Fourier 

(points) 

(averages) 

LAP 

Size 

1 

124.0 (133) 

129.5  (136) 

129.5(136) 

121 

Size 

2 

95.8 (108) 

120.3  (126) 

97.0(112) 

85 

Size 

3 

100.0(111) 

100.8  (113) 

112.8  (117) 

122 

Size 

4 

87.0 (95) 

82.3  (99) 

111.8  (117) 

105 

Table 

s  2b.  Average  values  of 

classification 

results  for 

single 

features . 


